home *** CD-ROM | disk | FTP | other *** search
/ SGI Developer Toolbox 6.1 / SGI Developer Toolbox 6.1 - Disc 4.iso / lib / mathlib / libblas / src_original / snrm2.f < prev    next >
Encoding:
Text File  |  1994-08-02  |  4.2 KB  |  154 lines

  1.       REAL             FUNCTION SNRM2( N, SX, INCX )
  2. *
  3. *     euclidean norm of the n-vector stored in sx() with storage
  4. *     increment incx .
  5. *     if    n .le. 0 return with result = 0.
  6. *     if n .ge. 1 then incx must be .ge. 1
  7. *
  8. *           c.l.lawson, 1978 jan 08
  9. *
  10. *     four phase method     using two built-in constants that are
  11. *     hopefully applicable to all machines.
  12. *         cutlo = maximum of  sqrt(u/eps)  over all known machines.
  13. *         cuthi = minimum of  sqrt(v)      over all known machines.
  14. *     where
  15. *         eps = smallest no. such that eps + 1. .gt. 1.
  16. *         u   = smallest positive no.   (underflow limit)
  17. *         v   = largest  no.            (overflow  limit)
  18. *
  19. *     brief outline of algorithm..
  20. *
  21. *     phase 1    scans zero components.
  22. *     move to phase 2 when a component is nonzero and .le. cutlo
  23. *     move to phase 3 when a component is .gt. cutlo
  24. *     move to phase 4 when a component is .ge. cuthi/m
  25. *     where m = n for x() real and m = 2*n for complex.
  26. *
  27. *     values for cutlo and cuthi..
  28. *     from the environmental parameters listed in the imsl converter
  29. *     document the limiting values are as follows..
  30. *     cutlo, s.p.   u/eps = 2**(-102) for  honeywell.  close seconds are
  31. *                   univac and dec at 2**(-103)
  32. *                   thus cutlo = 2**(-51) = 4.44089e-16
  33. *     cuthi, s.p.   v = 2**127 for univac, honeywell, and dec.
  34. *                   thus cuthi = 2**(63.5) = 1.30438e19
  35. *     cutlo, d.p.   u/eps = 2**(-67) for honeywell and dec.
  36. *                   thus cutlo = 2**(-33.5) = 8.23181d-11
  37. *     cuthi, d.p.   same as s.p.  cuthi = 1.30438d19
  38. *     data cutlo, cuthi / 8.232d-11,  1.304d19 /
  39. *     data cutlo, cuthi / 4.441e-16,  1.304e19 /
  40. *
  41. *     .. Scalar Arguments ..
  42.       INTEGER                          INCX, N
  43. *     ..
  44. *     .. Array Arguments ..
  45.       REAL                             SX( 1 )
  46. *     ..
  47. *     .. Local Scalars ..
  48.       INTEGER                          I, IX, J, NEXT, NN
  49.       REAL                             CUTHI, CUTLO, HITEST, ONE, SUM,
  50.      $                                 XMAX, ZERO
  51. *     ..
  52. *     .. Intrinsic Functions ..
  53.       INTRINSIC                        ABS, FLOAT, SQRT
  54. *     ..
  55. *     .. Data statements ..
  56.       DATA                             ZERO, ONE / 0.0E0, 1.0E0 /
  57.       DATA                             CUTLO, CUTHI / 4.441E-16,
  58.      $                                 1.304E19 /
  59. *     ..
  60. *     .. Executable Statements ..
  61. *
  62.       IF( N.GT.0 )
  63.      $   GO TO 10
  64.       SNRM2 = ZERO
  65.       GO TO 140
  66. *
  67.    10 ASSIGN 30 TO NEXT
  68.       SUM = ZERO
  69. *
  70. *        begin main loop
  71. *
  72.       IX = 1
  73.       IF( INCX.LT.0 )
  74.      $   IX = 1 - ( N-1 )*INCX
  75.       I = IX
  76.       NN = IX + ( N-1 )*INCX
  77.    20 GO TO NEXT( 30, 40, 70, 80 )
  78.    30 IF( ABS( SX( I ) ).GT.CUTLO )
  79.      $   GO TO 110
  80.       ASSIGN 40 TO NEXT
  81.       XMAX = ZERO
  82. *
  83. *        phase 1.  sum is zero
  84. *
  85.    40 IF( SX( I ).EQ.ZERO )
  86.      $   GO TO 130
  87.       IF( ABS( SX( I ) ).GT.CUTLO )
  88.      $   GO TO 110
  89. *
  90. *        prepare for phase 2.
  91. *
  92.       ASSIGN 70 TO NEXT
  93.       GO TO 60
  94. *
  95. *        prepare for phase 4.
  96. *
  97.    50 I = J
  98.       ASSIGN 80 TO NEXT
  99.       SUM = ( SUM / SX( I ) ) / SX( I )
  100.    60 XMAX = ABS( SX( I ) )
  101.       GO TO 90
  102. *
  103. *        phase 2.  sum is small.
  104. *                  scale to avoid destructive underflow.
  105. *
  106.    70 IF( ABS( SX( I ) ).GT.CUTLO )
  107.      $   GO TO 100
  108. *
  109. *        common code for phases 2 and 4.
  110. *        in phase 4 sum is large.  scale to avoid overflow.
  111. *
  112.    80 IF( ABS( SX( I ) ).LE.XMAX )
  113.      $   GO TO 90
  114.       SUM = ONE + SUM*( XMAX / SX( I ) )**2
  115.       XMAX = ABS( SX( I ) )
  116.       GO TO 130
  117. *
  118.    90 SUM = SUM + ( SX( I ) / XMAX )**2
  119.       GO TO 130
  120. *
  121. *        prepare for phase 3.
  122. *
  123.   100 SUM = ( SUM*XMAX )*XMAX
  124. *
  125. *        for real or d.p. set hitest = cuthi/n
  126. *        for complex      set hitest = cuthi/(2*n)
  127. *
  128.   110 HITEST = CUTHI / FLOAT( N )
  129. *
  130. *        phase 3.  sum is mid-range.  no scaling.
  131. *
  132.       DO 120 J = I, NN, INCX
  133.          IF( ABS( SX( J ) ).GE.HITEST )
  134.      $      GO TO 50
  135.          SUM = SUM + SX( J )**2
  136.   120 CONTINUE
  137.       SNRM2 = SQRT( SUM )
  138.       GO TO 140
  139. *
  140.   130 CONTINUE
  141.       IF( I.NE.NN ) THEN
  142.          I = I + INCX
  143.          GO TO 20
  144.       ENDIF
  145. *
  146. *        end of main loop.
  147. *
  148. *        compute square root and adjust for scaling.
  149. *
  150.       SNRM2 = XMAX*SQRT( SUM )
  151.   140 CONTINUE
  152.       RETURN
  153.       END
  154.